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ABSTRACT 


An analytical approach to the problem of replenishing task forces 
at sea 1s investigated, using a random walk model. The type of 
replenishment operations considered consists of two replenishment 
ships of different types, each replenishing L combatant ships of 
different types. The combatant ships are initially distributed into 
queues of length M and N between the two replenishment ships. As each 
combatant ship completes its first replenishment, it then enters the 
queue of the other replenishment ship. Thus each combatant ship is 
replenished by both replenishment ships. Replenishment times for one 
combatant ships are assumed to be independent, positive valued random 
variables whose distributions are dependent upon both replenishment 
ship type and combatant ship type. The replenishment of the task 
force is completed whenever 2L replenishments, two for each combatant 
ship, are completed. The total replenishment time is shown to be the 
maximum of certain partial sums of the individual replenishment times. 
A general expression for the distribution function of the total re- 
plenishment time is sought by various methods of analysis. A computer 
simulation is used to gain a better understanding of the form of the 
distribution function of total replenishment time for a specific under- 


way replenishment operation. Some examples are discussed. 
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I. THE UNDERWAY REPLENISHMENT OPERATION 


The underway replenishment operation is the means by which needed 
supplies are provided to Navy fleet units while underway at sea. It 
provides these units with the capability of operating at great distances 
from Navy shore facilities, enables them to remain on station for ex- 
tended periods, extends their period of effective combat readiness, 
and limits diversion from their assigned mission. The supplies pro- 
vided may include, but need not be limited to, petroleum products, 
missiles, ammunitions, stores, and provisions. However, no matter 
what is to be transferred between the transferring replenishment ships 
and the receiving combatant ships, several points are apparent, the 
underway replenishment operation increases the vulnerability of the 
engaged units, restricts their mobility, and diverts the involved 
combatant ships from their principal missions. The Navy is continually 
Striving for the development of techniques and methods which will 
reduce the limitations imposed on fleet operations by the necessity 
for underway replenishment. Principally, the Navy desires to reduce 
the vulnerability of the engaged units, increase the unit mobility 
while engaged, reduce the time required to complete replenishment, 
and limit interference with the assigned missions of the combatant 
units. 

To attain these goals many improvements of a material nature have 
been accomplished. The introduction of faster replenishment ships, 
such as the combat stores ship (AFS), reduces the vulnerability and 
increases the mobility of the engaged units by allowing underway 


replenishments to be conducted at higher speeds than before. The 


introduction of multi-commodity replenishment ships, such as the fast 
combat support ship (AOE) and the fleet replenishment tanker (AOR), 
reduces the total number of replenishments required since several needs 
can be filled during a single replenishment. The introduction of the 
helicopter as a transfer agent affords greater separation of task force 
units, reduces the demands for alongside replenishments, and minimizes 
interference with the assigned mission of the involved combatant ship. 
The introduction of new wire-highline transfer systems, such as the 
fast automatic schuttle transfer cargo handling system (FAST), provides 
a higher rate of cargo transfer and an increased load capacity, thereby 
reducing the required time alongside. In addition, internal cargo 
handling improvements in both the receiving and transferring ships 
have resulted in reduced breakout and strikedown times, thus again 
resulting in reduced replenishment times. All these improvements have 
had one element in common. Their accomplishment has resulted in a 
reduction in the actual replenishment time for an individual ship, and 
a subsequent reduction in the total replenishment time for the involved 
task force. 

Improvements have not been limited to those of a material nature. 
The underway replenishment operation has also been viewed as a problem 
in operations analysis. The objective of such an approach is a better 
understanding of the underway replenishment operation as viewed in an 
analytical manner. 

The underway replenishment operation can be viewed in many analytical 
contexts, but more readily lends itself to representation as a multi- 
Server queueing model with a finite number of customers. The number of 


servers can be established by the number of replenishment ships involved, 


10 


while the number of customers can be established by the number of ships 
in the task force scheduled for underway replenishment operations. It 
should be noted that, under normal conditions, a definitive sequence of 
individual replenishments has been established prior to the commencement 
of the underway replenishment operation. This is done as a method of 
optimizing replenishment ship utilization, while hopefully minimizing 
total replenishment time for the task force. Thus, while the underway 
replenishment operation may be viewed in terms of a queueing problem, 

it has definite job-shop scheduling problem overtones. 

One aspect of a typical queueing model must be ignored for the 
underway replenishment operation. It is normal to consider all service 
times to have come from the same probability distribution in a typical 
queueing model. Since many different types of replenishment and com- 
batant ships are involved in a typical underway replenishment operation, 
this seems to be an unreasonable assumption for any realistic model of 
the operation. Analysis of actual underway replenishment data has 
Supported this contention and has yielded approximations for the 
appropriate distribution functions for various typical underway 


replenishment operations. 
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II. THE UNDERWAY REPLENISHMENT MODEL 


The basic underway replenishment model that was developed by 
Waggoner [1] and later expanded by Patterson [2] is the domain of 
inquiry for this thesis. This model pen: that two replenishment 
ships and L combatant ships _ scheduled to conduct underway replenish- 
ment operations; a total of 2L replenishments are to be conducted, L by 
each replenishment ship. Prior to the commencement of the replenishment 
Operation, the L combatant ships are given a specific replenishment 
Sequence and assigned to queues of neck M and N, where M+ N = L, 
associated with the first and second replenishment ships, respectively. 
For each combatant ship this means that a definitive order of replenish- 
ment has been established and that each combatant ship has initially 
been assigned to one of the two queues. In addition, this means that 
the individual combatant ship both succeeds and precedes the same ships 
in each queue and that this queueing sequence is inalterable. 

The movement of the combatant ships during the underway replenishment 
operation is envisioned as a sequence of steps in a random walk. This 
random walk is performed on an Lx L grid on which the coordinates of 
each grid node represent a possible combination of the numbers of com- 
batant ships replenished by each replenishment ship, i.e., the node 
(i,j) on the grid indicates that i combatants have been replenished by 
the first replenishment ship and j by the second replenishment ship. 

It is obvious from the replenishment sequence that certain grid Rode: 
represent infeasible combinations; it is infeasible, for example, for 
the first replenishment ship to have completed M + k replenishments 


unless the second replenishment ship has completed at least k, where 
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k = 1,2,...,N. This infeasibility can be stated in constraint equations. 
If the abscissa, Zi» of the grid denotes the number of replenishments 
completed by the first replenishment ship, while the ordinate, 2,,of the 
grid denotes the number of replenishments completed by the second re- 


plenishment ship. The following equations constitute the necessary 


constraints: 
Zs = 0 
Z =M+ Zo 
Z, = L 
Z5 = 0 . 
Z. =N+ Z 
Z. = L 


All feasible nodes either lie on or within the boundaries established 
by the intersection of these constraint equations. Figure 1 depicts 
an Lx L grid with the constraint equations superimposed. 

In order to establish a common terminology for later discussion 
and investigation of the model, various definitions and theorems 
Stated by Patterson [2], are presented. 

Definition l 
.The CONSTRAINT BOUNDARY of the random walk grid is defined as follows 
(see Figure 1): 


a. the UPPER BOUNDARY consists of those nodes along the line 


Z. = L with abscissa k, where M < k < L-1; 

b. the RIGHT BOUNDARY consists of those nodes along the line 
Z, ob With OFC ENate ole whe he | Nicmiimces b= 4), 

c. the LOWER BOUNDARY consists of those nodes along the line 
Z) = M+ Zo» where O < 2, es ak setae 


3 


LEFT BOUNDARY UPPER BOUNDARY 


(L,L) 
(L,M) (it= leas) x 
(L,L-1) 





RIGHT 
BOUNDARY 





(s ,vtl) LOWER 
BOUNDARY 


Underway Replenishment Model 


Random Walk Grid 


FIGURE 1 
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d. the LEFT BOUNDARY consists of those nodes along the line 


Z, = N+ Z 


9 i. where 0< Z 


is este 

Definition 2 

A STEP is a horizontal or vertical line between adjacent nodes of the 
grid; a PATH is a sequence of non-negative steps which starts at the 
origin and terminates at (L,L) and does not cross the constraint 
boundary; and X. and we denote the times required to complete the ith 
horizontal and jth vertical step, Segreetiiee iby, of any path. 

From Definition 2 it is evident that Xx. is the ith replenishment 
time for the first replenishment ship, while me is the jth replenish- 
ment time for the second replenishment ship. S 
Definition 3 
The set of all feasible paths from the origin to (L,L) is decomposed 
into four mutually exclusive and exhaustive subsets: 

A: the set of all paths that pass through the node (L-1,L) and 
do not pass through any nodes on the lower boundary; 

B: the set of all paths that pass through the node (L-1,L) and 
pass through at least one node on the lower boundary; 

C: the set of all paths that pass through the node (L,L-1) and 
do not pass through any nodes on the left boundary; and 

D: the set of all paths that pass through the node (L,L-1) and 
pass through at least one node on the left boundary. 

In addition, subsets B and D can be further decomposed: 

By: the subset of all paths in B which last pass through the 
lower boundary at node (M+ k,k) for k = 0,1,...,N-1 

D,: the subset of all paths in D which last pass through the 


h 
left boundary at node (h, Nth) for h = 0,l1l,..., M-l. 
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Theorem l 
The total replenishment time is the actual time required to complete a 
path through the grid and is equal to one of the following values: 


a. the total length of time required for any path in A is 


b. the total length of time required for any path in B, is 


k 
' L k+1 
eg EOE Ne Wa nani: 
i=Mte+1 j=1 


c. the total length of time required for any path in C is 


_d. the total length of time required for any path in D, is 


eel L 
Oy aes tens NV pil, oem. 
i=l j=Nth+1 


Waggoner and Patterson both assumed that the X. and a were inde- 
pendent, random variables from the family of exponential probability 
distributions. Waggoner tacitly assumed that the replenishment rate 
was dependent only upon the replenishment ship type and that replenish- 
ment rates X} and wy could be used for all replenishments conducted by 
the first replenishment ship and the second replenishment ship, res- 
pectively. Patterson generalized this assumption to the extent that 
the replenishment rate was allowed to be dependent upon both the 
replenishment ship and the combatant ship types involved. Thus re- 
plenishment rates hs and MW, were assumed for the ith replenishment by 
the first replenishment ship and the jth replenishment by the second 
replenishment ship, respectively. It should be noted here, however, 


that Patterson developed Theorem 1 without referral to his assumption 
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of exponentiality. Thus Theorem 1 holds without any assumption being 
made regarding the distribution of the individual replenishment times. 
An interesting comparison can be made regarding the investigative 
approaches of Waggoner and Patterson. Waggoner approached the model 
by considering the first step taken and the subsequent sequence of 
steps of shortest duration taken during the random walk. Thus the 
exponentiality assumption, by means of the memoriless property, 
insured that the duration of each step continued to be exponentially 
distributed. This greatly simplified his derivation of the Laplace 
transform for the distribution function-of the total replenishment 
time. On the other hand, Patterson reversed the sequence and con- 
sidered the last step taken and the subsequent sequence of steps of 
longest duration in the reverse order. In this case, as noted above, 
the exponentiality assumption was not applied. Patterson, then using 
the exponentiality assumption, was able to’develop the Laplace trans- 
form for the distribution function of total replenishment time only 
for the case of two combatant ships initially distributed one to each 


replenishment ship. 


ly 








III. TOTAL REPLENISHMENT TIME 


The previous section has intimated that total replenishment time 
can be one of certain feasible values. It was, in fact, stated that 


: ; ; k . 
total replenishment time was either Ty» Oner (for some integer k 


B 


where 0< k < N-1), or T 


C opened k 


B (for some integer h where 0 < h < M-l). 


In this section, by means of the following theorem and proof, it is 
shown that total replenishment time is actually the maximum of those 
feasible values. 

Theorem 2 

The total replenishment time, i.e., the time required to move from node 
(0,0) to node (L,L) of the random walk grid, is the maximum of all 


feasible values. Let T denote total replenishment time, then 


k 


: for h = 0,1,...,M-1). 


T = max (T,3; T for 1 =.0el.. .-.o3N=13 Toa; ik : 


D 


Proof 


a. Assume that the path taken was from A. Thus it must be shown 


and T, eal a for 


that T, > T Gore fe 2 Ooi, .- Nelle ae Soa - 


B A ene 
eel... . M1. 
1) Because the path is in A, it must not touch the lower 


boundary. Thus it must reach the line Zo =k +1 (see Figure 1) before 


it reaches the line Z, = M+ k for any k O,1,...,N-l. This implies 


1 

that 

Mtk k+1 

hee 2 Vo for any k=O, 1 ..N-1. 

, i j 

i=l jal 

L 
Now by adding the partial sum o X. to both sides of the above 
i=Mtk+ 1 


inequality, it is seen that Ty > ay for any k = 0,1,...,N-l. 
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2) Because the path is in A, it must reach the upper boundary 


before it reaches the right boundary. This implies that 


Therefore, T, = To: 


3) Let the path in A first touch the line Z, = N+h (see 


2 
Beareure 1).at the node (r,Nth) where h < ten E-litorsany ti) 0st ioe. 
Since the line Z. = N+ h must be reached vertically, the (r + 1)th 
horizontal step must have begun before the node (r,Mth) is reached. 
The (Mht+-1)th vertical step is either just beginning (htl < r < L-1) 
at node (r,Mh) or will begin on completion of the (r + 1)th horizontal 
step ( r =h). Since the upper boundary must be reached before the 
right boundary, this implies that either 
L IL 
» X. = De LT ine or on 
i=rt1 j=Nth+1 J 
L L 
> X. = ye Ys tee hh: 
i=rt2 j=Ntht+1 J 
In either case, the above inequalities imply that 
L L 


oy X. > Due LOG (any enna O Ma 
i=h+2 j=Meht1 J 


h+1 
Now if the partial sum > X. is added to both sides of the last 
i=l 
inequality, it can be seen that T, = jie ror any n= Of bv 
b. Assume that the actual path taken was from Be where k is fixed 
and O< k= N-l. Thus it must be shown that co > TA3 Te = ae Ons 
oer ki Oe on =< Na) 8 : etietmeatid f - > on ‘ homens 0, ,. 22 Mole 
> ae pas b B CG? B D > > > 
1) Because the path is in Bes it must reach the line zy =Mt+k 
(see Figure 1) before reaching the Z, = ko ee Chus 
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j=l J is. ? 
L 
Now by adding the partial sum 2 X. to both sides of the above 
i=Mt+kt+ 1 


inequality, it can be seen that ie = Ty. 


2) Because the path is in B,, it must pass through the node 


es 
(Mtk,k). Thus the (Mtkt1)th horizontal step begins at the same time 


as the (kt2)th vertical step. Since the vath cannot pass through any 


node (Mtu,u) where k < u < N-1, it must reach the line Z. =u+tl 
(see Figure 1) before reaching the line zy = Mtu. These conditions 
imply that 
Mtu X > ut] 
s i Ds, HOG: eu ale 
i=Mtk+1 j=kt2 
L k+l 
Now by adding the two partial sums > xX, and & Y,. to both sides 
: L ‘ 
i=Mtut1 j=l 
of the above inequality, it can be seen that he = Tee for k<u< N-l. 
3) Let the path in BL first touch the line Z, = v+l (see 


Figure 1) at the node (s,v+l) where 0 < s < Mtv and 0<v<k. The 
(vt+2)th vertical step begins at this node, while the (s+l)th horizontal 
step either just begins (s = Mtv) or has already begun (0 < s < Mtv-1l). 


Since it is necessary that the line Z, = ktl (see Figure 1) not be 


Z 


reached until after the line Z) =M+k is reached, this implies that 


either 
k+1 Mtk 
ee. > xX. ifs=Me¢+v or 
jevt2 J isstl 
k+l Mt+k 
SS ye x X. Hye (Ce, ae ily 
javt2 J iss42 


In either case, this implies that 
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k+1 Mtk 
ee y X, for O< v<k. 
javt2 2 i=Mtvt1 


Jarl L 
Now if the partial sums % Y,. and > XX. are added to both sides of 
j=1 2 i=Mtk+1 


the last inequality, it can be seen that is Be es fom Ol ay aie 


4) Because the path is in B the (k#2)th vertical step and 


ic? 
the (Mtk+1)th horizontal step begin at the same time. Since the 


upper boundary must be reached before the right boundary, this implies 


that 
L L 
~ xX, > awe a 
k=Mtkt1 j=kt2 J 
k+1 
Now by adding the partial sum % Y. to both sides of the above 
jel ? 
inequality, it can be seen that ts > To: 
5) Because the path is in By it cannot reach the line Z, = k+l 
(see Figure 1) before reaching the line zy = Mtk, which implies that 
k+1 Mtk 
Die SS X. 
j=l J isl 
This in turn implies that 
k+1 pall 
(1) yy Re 2 ek torngany one 10s ML. 
: j it 
j=l i=1 


Further, let the path yee) BL first reach the line Zo =N+h (see 
Figure 1) for any h = 0,1,...,M-1 at some node (w,Mth) where Mtk < w < L. 


Since the path reaches Z, = N+ h through a vertical step, the (Mhtl)th 


Z 
step is just beginning at node (w,NMth) while the (wtl)th horizontal step 
has already begun. Since the upper boundary must be reached before the 
right boundary, it is necessary that 

L L 


u X, > De een eae eco aemler 
i=wtl j=N+h+1 
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This in turn implies that 
L L 
(2) 2 Xx. > 2 ceeok sckOl gamiy aig => Ole ees holes 
i=Mtkt1 j=Meht1 J 
Now if the inequalities (1) and (2) above are added together, it can be 
k h 
seen that Ta > Ty Forrhm— On. Lo 22 2M—1. 


Ge) lif the path takenswas in ©, Po tcelancer thane eicmother 


C 
feasible values. This is proven analogously to part a. 


d. If the path taken was in Di» for some h =205 1... Mit Th 


is larger than the other feasible values. This is proven analogously 
Go part b. 
This concludes the proof of Theorem 2. 

By way of illustrating the result of Theorem 2, an example where 
constant individual replenishment times are involved is discussed. 
It is assumed that x. = A and ve = B for ig). = 162,42. soo aaene 
general case of L combatant ships, distributed into queues of length 
M and N, associated with the first and second replenishment ships, 
respectively. If neither N = O nor M = OQ, the following are the 


feasible values for T: 


T, = IA 
T, = (N-k)A + (kct1)B for k = 0,1,...,N-1 
T. = 1B 

et 5 A, 2 OR Hor 0, ae 


T, the maximum of all feasible values, is clearly equal to either LA 
or LB, depending upon the values of A and B. If A >B, then T = IA;. 
if B >A, then T = LB. If N = 0, all feasible paths are in D and the 
following are the feasible values of T: 


T = (EDV EO ac hh S Onl oie 


af 








T is clearly equal to either (A + LB) or (IA + B), depending upon the 
values of A and B. If A>B, then T = (IA +B); if B >A, then 
T = (A+ LB). If M = 0, the preceding argument holds if the set D 


is replaced by the set B. 
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IV. IDENTIFICATION OF THE GENERAL PROBLEM 


In the preceding section it has been shown that the actual total 
replenishment time, denoted T, is the maximum of all feasible values. 
Thus, while the exact composition of T as expressed in terms of the 
individual replenishment times is known, the general form of its dis- 
tribution function as expressed in terms.of those for the individual 
replenishment times remains unidentified. In regard to this last 
Statement, in this or any succeeding section the following notation 
has been used: the distribution function of X. is denoted F., while 
that of the ue is denoted ae if densities are assumed, they will be 
denoted by the appropriate small letter with subscript. 

The composition of T and the distribution function of T when Ga Sayan 
while of little practical interest, are presented below as an introduction 
to the analytical techniques used in the attempt to identify the general 
distribution function of T, when L > 3: 


a. if L= 1 with M = 1 and N = 0, then T = Xy sF vy and its dis- 


tribution function is the convolution of Fi and Gi; 


b. if L = 2 with M = 2 and N = 0, then T = X, + max (X, .Y,) + Y 


1 Z. 


and its distribution function is the convolution of F and G 


Pe Ho. Ze 


where Ho denotes the distribution function of max (X,Y), Le) 6 


5&3 


ec. if L =2 with M = 1 and N = 1, then T = max (X,5¥,) + max (X,Y) 


and its distribution function is the convolution of Hy and Ho 9» where 


H.; denotes the distribution function of max (X,Y. ), 5 A FG, Lot 


ae yp 
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It should be noted that for all cases investigated, either above 
or beyond this point, an assumption that M > N has been made. The 
symmetrical nature of the model has insured that any results obtained 
here are valid upon the interchange of M with N, X. with ee and F. 
with a for the case when M < N. 

The general distribution function of T can be represented by the 
following probability statements: 


a.8 Le N = 70>) Ehen 


+1 L 
Ge <<) (yeas Pye ee hoe (n= Oi . +1); 
i=l jehbt1 J 
ey rhs yh ee il 5 tela 
Jie IT k+1 
P@ <t)=P( 0X, <t; eee) 2 eee Sains. leds 
i=l i=Mtkt1 j=l 
L h+1 L 
DS ope: Comer SO ee iN = pil, .M-1). 
j=1 J i=l jaNthtl 


Viewing the summation terms involved in the two probability state- 
ments as random variables, denoted as To: it can be noted that when 
N = 0 L such random variables are present, while if N>1L+ 2 such 
random variables are present. Since it was desired to use the technique 
of transformation of variables to simplify these probability statements, 
it was obvious that additional random variables had to be introduced in 
order that the one for one transformation demanded by the technique 
could be attained. Since the transformation was to be from the Xx. and 
the me into the To» L artificial random variables would have to be 
introduced if N= 0, while L - 2 would suffice if N>1. These 
artificial random variables were selected to fulfill two purposes: 


primarily, they insure a one-to-one transformation, while secondarily, 


Zo 








they have no effect on the appropriate probability statements. For 


L 
the case of N = 0, random variables of the form }%& X. for m= Le2sre 
1=m L 
are introduced. These summations are specifically less than X. et Yh: 
L i=] 
For the case of N > 1 random variables of the form © X, for 
L i=m * 
m= 2,3,..., M and »} Y. for n = 2,3,...,N are introducedtatnece 
jen J L it 
summations are specifically less than % X. and }%Y Y. respectively. 
i=l j=l 


Thus for either case the introduction of these particular summations 
into the appropriate probability statement have no effect. 
The whole set of new random variables is defined below: 


a. if N=0, 


fe 
T= 3 Xe) form = 1,2, »L 
i=m 
h+1 L 
z eee - >) Bee ror be —70 5 le bale 
Ith+ 1 1 k see 
b TERN eo .. 
ib 
T= 2 Meeeor me = 16 2, M 
m : i 
i=m 
ic k+1 
Lt = Kk. ue VerotiK — 0) banal 
bs ae iets) adel 
L 
Then = ZY. for n = 1,2, ,N 
Jao 
sae Jb 


TeNth+1 = je Xx. + 2 Ysfor-h = 041.2 el 


The Jacobian matrices of these transformations are shown in Figures 2, 

3 and 4. The value of their determinants is the subject of Lemma 1. 
Lemma 1 

The Jacobian determinant of the transformation matrices shown in Figures 


2, 3 and 4 is equal to 1, unless L = 2 with M = 1 and N= 1, i.e., if 


26 








2L 


JACOBIAN MATRIX OF T TRANSFORMATIONS WHEN N= O 


1 2 L-1 L 1 Z L-1 
1 1 1 1 0 0 0 
0 1 1 1 0 0 0 
0 0 0 1 0 0 0 
1 0 0 0 1 1 1 
1 L 1 1 0 0 0 
FIGURE 2 


27 





JACOBIAN MATRIX OF T TRANSFORMATIONS WHEN N = 1 


1 2 M-1 L 1 2 
Ty Ie 1 1 1 1 0 0 
T, 0 1 1 1 1 0 0 
Ty 0 0 0 1 I 0 0 
qT, 0 0 0 0 it 1 0 
Tiel 0 0 0 0 0 1 1 
T 42 It 0 0 0 0 0 i 
Tox. 1 1 I 1 0 0 0 

FIGURE 3 
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| 5| denotes the determinant of the Jacobian transformation matrix J, 
then [s| = 1. 

Proof 

a. Rows 1 through L +N have all O's to the left of the main 
diagonal elements and the main diagonal elements are equal to l. 

b. Rows L+N+ 1 through 2L-1 can be reduced to similar form by 
the application of the formula 


' = - eS = 
RY; aa eRe a wg go> Ba F Brag com ten cee 


where R. denotes the ith row before application of the formula, and 
R', afterwards. 
c. Row 2L can be reduced to a similar form by application of the 
formula 
R' or = Ror - Ry if N = 0, 
By ot te Ge ee ae Bee Me NS, x 


i] = = = : 
Te Sere Teer" are een erp i eae 


The case for L = 2 with M = 1 and N 


l was specifically excluded 
from Lemma 1 since the determinant of its Jacobian transformation 
matrix equals zero. This results from the inherent dependence of the 
involved summations. It is noted that this is the only case in which 
no artificial random variables can be introduced. Of course, the 
distribution function of T in this case was given earlier in this 
Section. | 
A second transformation of variables was made at this point to gain 

further simplification. The following random variables were defined: 


Vie for m= 1,2,...,2L-1 
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The determinant of the resulting Jacobian transformation matrix equals 
1 since the matrix is triangular with all main diagonal elements equal 
to 1, i.e., all elements below the main diagonal are 0's. 

At this point it was decided to assume that the distribution 
functions of the individual replenishment times were continuous and 
possessed non-vanishing density functions when the Xx. or Me are ee. 


than 0. The joint density function of the Xx. and the = can thus be 


written 
L L 
£ <x.) Tee.) x. S.0 andiy es =a0 
Le See 2 aw od 
£(X) 5+ + XpsYyo+- + Vz) = for 1,j=l,2.6een 
0 otherwise, 


since it was previously assumed that the xX. and the ue were independent 
random variables. 
Since the appropriate Jacobian determinants are equal to l, the 


joint density function of the We can also be written in the same 


manner 
L L 
f.(x.) 1 g.(y.) x. > 0 andy. > 0 
i=l * * j=l Jod ‘ J 
B(Vy5Vo9+++2Voq) = 9) Salas beets ad I ae © 
0 otherwise. 


The region in V-space for which gv, Voz) does not vanish can be 


sVooree 
represented in terms of the ve if the transformations relating the x. 


and the Me to the Ve are used. These are given below: 


a. if N= 0, 
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7 = el oF vee fOr] oe lee laa 

Di 
eo eee 

i1=L 

2L-l 
tr = - 2 Ve 

i=1 

b. if N=, 
X. = V Foc 1 = 12a Ma 
ee Same Viens, 
i, 
ge 
IF l L i=U+2 i 
iS = ite + Via for “j= 442, ,L-2 
ee Mh ea ees Sen 
ea 
ty, a 2V, - 2 vi ~ Wen oP oye 
{=2 
eee Ne 2, 
X. = V for i-= 12 
ie) dye) Vaca BOF be SI Opilsous Me” 
Se a al oP ea 
Dire 
ose yy Pay. 
. ; i=L+Mt1 * 

Y, 2 Vies fon 4G = 120m N21 
Yueh = Nieagh oy Vint fora — 0,1. »M-2 
ee eae Hope 

M 
“Lee Oe ee 


The above statements 


thus allow the original probability state- 


ments for the distribution function of T to be represented in terms 
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in the following form: 


ANN 
PCL < “t= PC EV, < t for m= 1,2,...,2L) 
i=m 


of the V 
m 


provided that the joint density function for the Va does not vanish 


As an example, the joint density function of the Va expressed wholly 


in terms of the V_, is presented below for the case N > 2 


M-l L-2 
ai £.W,) nf. , sty fy 70 stv, it 1 Vaan 
i=l] i=M 

IRE N-1 L-2 


ff. (-v,tv_t+ u Vs) 1 8 ) Tt Save +v_..) 
glee i=LtNeL st Lj jen j° jtl-N Ltj 


B(V)5Vo0-++2Vo7) = 
Bp OM i Van- eto By cian 


Ny condition (1) below exists 


0 otherwise. 


Condition (1) is expressed below: 


Ve, ea) ye reren y Eoall 
V Varied OO 1 = M, 2.25 L-2 
oy eS 10 


It is evident that the above discussion has not concluded in a 


truly desirable result, namely a general form of the distribution 
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function of T that would facilitate ease of computation. The resulting 
forms, those containing the Th and the har do, however, exhibit prop- 
erties that may prove useful in future attempts in this area. As the 
matter now stands, not even a relatively simple problem, such as that 
when L = 3 with M = 2 and N = 1, can be directly resolved. Appendix A 
depicts this particular problem as an illustration of the general dif- 
ficulties encountered. In retrospect, it can be stated that the 
dependence relationships between the T > or those between the Rene 
impose restrictions that cannot be overcome without Some convenient 
manner of handling the joint distribution function of dependent random 


variables, such as that developed in this section. 
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V. A SIMULATION OF THE UNDERWAY REPLENISHMENT OPERATION 


The preceding section, in which analytical methods were used, 
failed to produce the distribution function of T in a form feasible 
for specific computations. It was therefore decided to study the dis- 
tribution of T by computer generation of T realizations for a specific 
case of underway replenishment. The primary intent of such a study was 
to gain a — understanding of the form of the distribution function 
of T for typical underway replenishment operations. 

The general method used made direct use of the results of Theorem 
2, i.e., that T is the maximum of all feasible values. The underway 
replenishment operation was not directly simulated. Instead, the 
individual replenishment times were simulated and the resultant values 
of T calculated, using the results of Theorem 2. 

For the specific underway replenishment operation, it was decided 
to use an attack carrier task force replenishment operation. The task 
force was selected to be composed of eight combatant ships, an attack 
@ircraft carrier (CVA), a light guided missile cruiser (CIG), two 
guided missile frigates (DIG), and four destroyers (DD). The two ship 
replenishment group was selected to be composed of a fleet oiler (AO) 
and an ammunition ship (AE). This particular selection was made because 
approximations for the distribution functions of the individual replenish- 
ment times for the above type ships were presented by Besecker in [3]. 
Thus Erlang distribution functions were assumed for the individual 
replenishment times; the specific gamma parameters, denoted k and}, 
are shown in Table I. The general form of the Erlang distribution 


function is given in Appendix B. 
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PARAMETERS FOR ERLANG DISTRIBUTION FUNCTIONS 
FOR 


INDIVIDUAL REPLENISHMENT TIMES 


REPLENISHMENT SHIPS 





C 

O CVA ll O09 3 0.03 
MS 

B oH CLG i 0.11 2 0.02 
A I 

I ?P DLG 7 Orel 2 0.02 
A Ss 

N DD 4 0.05 3 Oe tels 
dk 

TABLE I 


In accordance with the model, the AO and AE were designated the 
first and second replenishment ships, respectively. The AO is then 
associated with the Xs; While the AE is associated with the ae It 
was decided to establish a constant replenishment sequence at the AO. 
This sequence is CVA, DLG, DD, DD, CLG, DLG, DD, and DD. The replenish- 
ment sequence at the AE would depend upon the M and N combination being 
used for the generation of a T realization. The feasible sequences 
are shown in Table II. 

The individual replenishment times were simulated by the inverse 
transformation method discussed in Appendix B. The property that an 
Erlang random variable is the sum of k random variables, each expon- 


entially distributed with parameter 4, facilitates this method. To 
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insure that the random variates generated were from the assigned dis- 
tributions, chi-square goodness of fit tests were made at the 0.01 
level of significance. This method is discussed in Appendix B. For 
all 16 distributions of the individual replenishment times, the null 
hypotheses that the generated random variates came from the appropriate 


distributions were of course accepted. 





REPLENISHMENT SEQUENCE AT AE 






AHSAMS TM Sse zw 
feat (Gp) acd eal (om vem yen ep. 


TABLE II 


A special feature was provided in the implementation of the com- 
puter program to enable study of the dependence of T on the M and N 
combination being used. The program provided that upon the nth gen- 
eration of a T realization for a particular M and N combination, all 
ships would have the same individual replenishment times at both re- 
plenishment ships as for any other M and N combination, i.e., on the 


50th generation of a T realization for M = 3 and N = 5, the CLG would 
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have the same replenishment times at both the AO and AE as for the 50th 
generation of a T realization for any other M and N combination of the 
combatant ships. 

Five hundred T realizations were generated for each of the nine 
feasible M and N combinations. The random walk for a typical T 
realization is shown in Figure 5 as a random walk within the constraint 
boundary of the L x L grid and in Figure 6 as a random walk in two- 
dimensional--time space. The T realization obtained by performing the 
random walk agrees with that obtained through an application of 
Theorem 2. 

For each M and N combination, the generated data were used to con- 
struct frequency histograms of T realizations. Figure 7 is a typical 
frequency histogram. From these histograms, one basic observation 
could be made: the frequency histograms showed definite positive 
Skewness, i.e., they were skewed to the right. The obvious contention 
at this point was that the distribution function had the same character- 
istic shape as did those of the individual replenishment times, i.e., 
the shape of a gamma distribution function. 

A second observation could be made once the T sample means for the 
various M and N combinations were calculated. The results are similar 
to those cited by Milch and Waggoner in [4] for the expected values of 
T when the X. and He are exponentially distributed with parameters i 
and u, respectively. The minimal T sample mean occurred when M and N 
were equal and a range of M and N combinations exists in which there 
is only a slight variation of the T sample mean values. The T sample 
means are plotted against M and N combinations in Figure 8. 

It was decided to investigate the contention that the distribution 


function of T was from the family of gamma distributions. Estimates of 
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Individual Replenishment Times 





X = 123.35 yy = 150.71 
X, = 50.14 Y, = 88.46 
X, = 81.12 Y3 =P Zo, 
Xy, = 146.83 Yy, = 257.16 
Xe = 34.64 Y5 = 282.34 
Xe = ot. 9 Ye = 48.44 
Xo = 21.43 Y5 = 8.52 
Xe = 124.14 Ye es) Re 4 
T = 951.56 


A Typical Random Walk for the 
Specified Underway Replenishment Operation when M = 5 and N = 3 


FIGURE 5 
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Total Replenishment Time by AE (minutes) 


o (Completion at AE 
| at 951.56) 
O 





900 
| 
O 
800 
700 
600 
500 serra * (Completion at AO 


at 673.24) 





O 100 200 300 400 500 600 700 800 900 1000 


Total Replenishment Time by AO (minutes) 
xX: a completion of an individual replenishment 
by AO 
Oo: a completion of an individual replenishment 
by AE 
Note 
Data plotted is the same as used for Figure 5 
Total Replenishment Time by AO vs 


Total Replenishment Time by AE when M = 5 and N = 3 


FIGURE 6 
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Number of Observations in Class 


140 


130 


a 


5 fl) 521 17.61 220) 
we eae 11.01 15.41 19.81 


Total Replenishment Time (hours) 
Frequency Histogram of Typical T Observations 


when M = 3 and N = 5 


FIGURE 7 
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T Sample Mean Values (hours) 





M and N Combinations 


Relationship between T Sample Mean and Feasible 
M and N Combinations 


FIGURE 8 


42 





necessary gamma distribution parameters k and 4, denoted k andl were 
obtained by the method of moments as shown in Appendix B. The basic 
null hypothesis to be tested would be that the T realizations for a 
particular M and N combination came from a gamma distribution with 
parameters equal to k and Ne The chi-square method was choosen for 
testing goodness of fit. Despite the generally acceptable results, 
eight out of nine null hypotheses being acceptable at the 0.01 level 
OE significance, it was evident that better results could be achieved: 
By studying superimposed frequency histograms of the observed T 
realizations and the expected T realizations under the gamma distri- 
bution of the null hypotheses for both the accepted and rejected cases, 
two ovservations could be made: 

a. the peak ordinate value for the number of observed T realizations 
was higher than that of the peak ordinate value for the number of 
expected T realizations for all M and N combinations; and 

b. the mode, the absissca for the peak ordinate value, for the 
number of observed T realizations either occurred to the right of or 
approximately coincided with the mode for the number of expected T 
realizations for all M and N combinations. 

A typical superimposed frequency histogram is shown in Figure 9. 

From these observations, the following generalization could be 
made: 

the true density function of T has a maximum value that occurs 
above and to the right of a gamma density function with parameters 
k and A. 

Thus it was decided to use a method of estimating gamma parameters 


that utilized the maximum value of the gamma density function. This 
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Number of T Observations or Expected T Valuesin Class 





DG = 469 9.38 14.07 leo a a2Sen4 
De 7.03 11.72 16.41 21.10 


Total Replenishment Time (hours) 


Solid line: wtumber of T observations 
Dashed line: number of expected T values 


Note: Number of expected T values calculated using 
gamma distribution with parameters estimated 
by method of moments 


Superimposed Frequency Histograms 
of Observed and Expected T Values when M = 6 and N = 2 


FIGURE 9 
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method, termed the peak value method, is developed in Appendix B. 
Estimates of k and \, as determined by this method, are denoted k* and 
A*, respectively. Success was readily achieved by this method. The 
basic null hypothesis, that T for a particular M and N combination was 
distributed by a gamma distribution with parameters k*¥ and \¥*, was 
acceptable at the 0.01 level of significance by the chi-square method 
for all nine M and N combinations. The results, including those 
obtained by the methad of moments, are shown in Table III. It could 
be concluded that the distribution function of T could be successfully 


approximated by a gamma distribution function with parameters k* and \*. 


RESULTS OF METHOD OF MOMENTS AND PEAK VALUE ESTIMATION 





8/0 # rade 4.24 40.58 3.09 46.19 3 252 
7/1 12.28 4.08 36.95 3.01 ioe 3736 
6/2 11.83 3.76 BL ae SIS 44.63 a7 
cys M52 Bae CwiG 2 ee 48.66 LS 
4/4 ‘J 11.49 sly ‘Gy 3.63 47.70 5 
3/5 11.50 30 41.46 3.60 46.34 4.03 
2/6 11.53 B22 41.28 3.58 46.60 4.04 
Ly 11.95 3.97 35.91 3.00 45.06 oe 
0/8 146 Lian Sona3 2.86 46.93 3.49 


#: M and N combination that resulted in rejection of null 
hypothesis under method of moments testing 


TAB EES eet 
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APPENDIX A 


ILLUSTRATION OF PROBLEM WHEN L = 3 WITH M = 2 AND N= 1 
The elements of this particular problem are presented below in the 
context of Section IV: 
a. Considering the feasible paths, this problem involves random 


variables of the form: 


qT, = xy Br X, a X., 
ti = Kn oe Y 
B 3 1 
To = Yy A Y5 2 Y4 
ro = X a etc 
D 1 2 3 
ti =X, + X = ae 4 
D i 2 3 


The distribution function of T can be written 


P(T < t) = P(X +KjHX, St, XY, <t, YAY +¥, < t, 


nese, @ Ey MG ANGs ce fe). 


b. After the addition of the appropriate artificial random variable, 
the problem involves random variables of the form: 


ah xX age ne aoe 


1 1 2 3 


HT 


T, = Dian 58 


2 2 3 

T, = X., 2 Y 

qT), = Yy za) + Y, 
Ts = X, = Y, fe Y4 
Te. = Xy ats X, 25 Ys: 


46 





The distribution function of f can Desewrilteen 


Pa =e P(T, <befor aise)? war) 


c. After the transformation from the TO to the ee the problem 


involves random variables of the form: 


V,=%, 
Oo a 2 ee 

iy “Y, -¥, 
V, =X, +Y, 

Veo=  -X, +Y, 

Ve = X,tX, 4Y,. 


The distribution function of T can be written 
6 
PGi t) =P & Va CNPOr sig ly a One 
m=i 


d. Considering the inverse transformations from the Nos to the Xx. 


and the ae the problem involves the forms: 


Xi =. Vi 

X, = vy as V5 = wi 

X, = 4 a V3 te Ve ae Ve 
Y = vi 3 U7. 

Yo = Vi 1 Vy ot Win Ve 

Y3 =-2V, = V5 - My. oh Ve: 


e. The joint density function of the Va as wholly expressed in 


terms of the Mes can be written 
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f£ ea )f, (v ity tv, df 


3 ¢-V 1tV3tVctVe) ake 


8 (Vy V5 2V 4 2Vz V5 Vp) = 8, (v tv, )8, Vv tv, +v i tV 58s (- -2v, 
if condition (1) below holds 
0 otherwise 
Condition (1) is as follows 
vy > 0 
Vy a Vo a Mi, =, 
V4 + Ve apn = vs 
vy + Vy, 0) 
vy + Vo oF 7 oF Ve > 0 
V > 2V, + Vy ae 
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APPENDIX B 


RELEVANT THEORETICAL AND STATISTICAL FORMS 
1. The Gamma Distribution 
The gamma distribution function is given by 


x 


,k ot a AX dx 
F(x) = | 5“ p_S——* for x > 0 and k > 0, 


ae 
where I'(k) = [x e *~ dx. Let Us aia ela! po denote the mean, the mode, 


O 
and the variance of the gamma distribution, then 


u = k/\; 
m = (k-1)/A3; and 
s = Ne 


If k is a positive integer, then the resulting gamma distribution is 
also termed an Erlang distribution. The Erlang distribution can be 
shown to be that for the sum of k identically distributed independent 
random variables from an exponential distribution with parameter A. 
If k is of the form n/2 with n a positive integer and i} equals 1/2, 
then the resulting gamma distribution is also termed a chi-square 


arstribution. 


2. Determination of Erlang Random Variables by the Inverse Trans- 
formation Method 


Let F(x) denote the distribution function for a random variable X 
and let r denote a random number from the uniform distribution on the 
- -1 
interval 0 to l. If r = F(x), then x = F(r) a where F(r) denotes 


al 


n 
the inverse transformation of F(x). If X. is exponentially distributed 


* 
The transformation of a random variable X into a uniformly distributed 
random variable r = F(x) is a well known result in probability theory. 
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with parameter i, then 


= j- i: I = 
HGS) i) Mee 5 or if F(x.) =41r,, 
“AX. 
then r= l-e 1 ee 
thus NSH =] - tr. 5 bute it r. is a realization from the uniform dis- 


tribution on the interval 0 to 1, then l - Ee r,! is also la Trealizatron 
on the same interval. Thus this may be rewritten 


“XxX, ‘ ; 
e ier 3 solving now for Xs 


x .o= -" 1/\ olnere. 
i at 


k | 
nx = . 2, Xs then X has an Erlang distribution with parameters k and 


i=] 
4; thus 
k 
== 1/) Sein ri'3 which may be rewritten 
i=] 
k 
) aan ve in( mr,'), where 
aL 
i=] 
ml k k 
EG) = - l/r In(mr,") with c=) Woe 
3 1 : L 
i=] i=] 


3. Test for Goodness of Fit by Chi-square Criteria 

Let x. be. a realization from the sampled population; let 0, denote 
Bie number of realizations that fall into class j for j = 1,2,...,n, 
where n is equal to some arbitrary number of subdivisions of the range 
of the assumed distribution of Xs let Be denote the expected number 
of realizations in class j under the assumed distribution; then the 


Statistic 


n (E. - a 
——— S ie! ns Sees 
i=l BG 


has an approximate chi-squared distribution with parameter (n-m-1)/2, 


where m denotes the number of assumed distribution parameters estimated 
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with sample statistics. n-m-1l is called the degree of freedom associ- 
ated with X; for a choosen significance level the value of X can be 
compared with the appropriate table value such as found in Appendix 4 
of [5]; if X is greater than the table value, then the null hypotheses 
that the x, came from the assumed distribution is rejected; it is 
accepted otherwise. 
4. Estimation of Gamma Parameters by the Method of Moments 

Let T be a random variable from a gamma distribution with parameters 
k and A. Let op denote the expected value of T, T denote the T sample 
denote the variance of T, and 5? 


Z 
mean value, o denote the T sample 


Ak le 
n 
variance as calculated by l/(n-l) (,-T)’, where n is the number of 
i=l 
T sample values. Since 
Un = k/X, and 
Ze Z 
On = K/h”, 


A A 
it can be shown that k and \, estimates of k and 4, respectively, can 


be calculated by 


a> 
I 
a 
=a 
” 
NO 
Wy) 
~ 
Qu 


>> 
! 
I 
~ 
ee) 


2 , 
if the mean UW and variance Op are respectively estimated by the sample 


mean and variance. 
5. Estimation of Gamma Parameters by the Peak Value Method 

let T be a random variable from a gamma distribution with parameters 
k and X. Let Ur denote the expected value of it. 4 denote the T sample 
mean value, p denote the maximum value, i.e., the peak, of the relevant 
gamma density function, and m denote the mode of the distribution of T, 


A 
i.e., the T value at which p occurs. Let Lan De k* and A\* denote 
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estimates of Ups Ps k and X, respectively. Since p occurs at m and 


since m = (k-1)/\, p can be evaluated as 
— eerykt gtk 
P T'(k) 


Now by applying an approximation formula for [(k) as found in [6], 
i.e., 
C(k) = ./2n oe eid for k large and not an integer, 


p can be approximated as 


k-1 
a ee 
211k 


Since for large k, Gleyiey can be approximated by oie p can be 


approximated as 
A 


J/2Tik 





(1) Py 


A 
To estimate k and i, Lp and D will be used. Since Up = k/X and since 
p can be approximated by formula (1), by combining these equations and 


solving for k, it is easily seen that 

k~ 2np Hp. 
Now, knowing an approximation for k, an approximation of i can be 
obtained from A = K/Up- Now, by letting u = T and D = max ow, 
where max v is the number of T realizations in the modal observation 
class of the T sample frequency histogram, n is the number of T 


realizations in the sample, and w is the width of the modal class in 


the frequency histogram, k* and \* can then be evaluated as 


k* = ann T 
h® = k*/T. 


sys 
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An analytical approach to the problem of replenishing task forces 
at sea is investigated, using a random walk model. The type of 
replenishment operations considered consists of two replenishment 
ships of different types, each replenishing L combatant ships of 
different types. The combatant ships are initially distributed into 
queues of length M and N between the two replenishment ships. As each 
combatant ship completes its first replenishment, it then enters the 
queue of the other replenishment ship. Thus each combatant ship is 
replenished by both replenishment ships. Replenishment times for the 
combatant ships are assumed to be independent, positive valued random 
variables whose distributions are dependent upon both replenishment 
ship type and combatant ship type. The replenishment of the task force 
is completed whenever 2L replenishments, two for each combatant ship, 
are completed. The total replenishment time is shown to be the maximum 
of certain partial sums of the individual replenishment times. A 
general expression for the distribution function of the total replenish- 
ment time is sought by various methods of analysis. A computer simulation 
is used to gain a better understanding of the form of the distribution 
function of total replenishment time for a specific underway replenish- 
ment operation. Some examples are discussed. 
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